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Abstract 

This paper concerns kinematic helical dynamos in a spherical fluid body sur- 
rounded by an insulator. In particular, we examine their behaviour in the regime 
of large magnetic Reynolds number R„^, for which dynamo action is usually 
concentrated upon a simple resonant stream-surface. The dynamo eigensolu- 
tions are computed numerically for two representative single-roll flows using a 
compact spherical harmonic decomposition and fourth-order finite-differences 
in radius. These solutions are then compared with the growth rates and eigen- 
functions of the Gilbert and Ponty (2000) large Rm asymptotic theory. We find 
good agreement between the growth rates when Rm > 10'', and between the 
eigenfunctions when Rm > 10^. 



1. Introduction 

We consider a class of kinematic dynamos in which the magnetic field B 
is generated by the steady helical motion v of an incompressible, electrically- 
conducting fluid. Helical flows constitute some of the simplest and most efficient 
mechanisms for the excitation of a seed magnetic fleld, as witnessed in numerical 
simulation and exploited by laboratory experiments (Gailitis et al. 1987, Dudley 
& James 1989, Forest et al. 2002, Moss 2008, for example). Such flows are also 
widespread in astrophysical fiuid bodies, such as jets, the convection zones of 
stars and, possibly, liquid planetary cores, where they might appear as a slow 
meridional circulation (Giles et al. 1997, Gough and Mclntyre 1998, Olson et 
al. 1999, Haber et al. 2002, Hartigan et al. 2005). 

By far the simplest helical fiow is the Ponomarenko dynamo (Ponamerenko 
1973), which consists of a solid electrically-conducting cylinder of finite radius 
rotating with constant angular velocity fl and translating with a constant ax- 
ial velocity V, while its rigid (electrically-conducting) exterior remains at rest. 
The dynamo loop in this case comprises (a) generation of azimuthal magnetic 
field from radial magnetic field by the discontinuity in the rotation (the 'omega 
effect'), and (b) reciprocal generation of radial fleld from azimuthal fleld by 
magnetic diffusion. Note that the additional longitudinal shearing component 
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V is cnicial, as it draws apart oppositely directed field lines and prevents flux 
expulsion. Ponomarenko (1973) determined the magnetic field in this case ana- 
lytically, showing that the field is concentrated at the velocity discontinuity on 
the cylinder boundary. Helical dynamos with non-uniform V and and arbi- 
trary cross-section, such as in a conducting fluid, are naturally more complicated 
(but see the exact steady solutions of Lortz 1968, and Chen & Milovich 1984). 
Nevertheless, they share the same dynamo ingredients: differential rotation and 
magnetic diffusion. These ingredients are especially conspicuous at very large 
magnetic Reynolds number Rj^, a regime which can be probed successfully with 
asymptotic theory (Ruzmaikin ct al. 1988, Gilbert 1988, Gilbert & Ponty 2000, 
hereafter GP). In such models the importance of diffusion, in particular, is ev- 
ident from the localisation of potentially growing dynamo modes upon their 
critical surfaces, where the modes arc nearly convectcd by the flow. It is upon 
these surfaces that the modes exhibit the small-scale variation necessary for 
magnetic diffusion to work effectively and replenish radial magnetic field. Not 
every stream-surface can support a mode, but a 'resonance condition' selects 
the critical surface (or surfaces) upon which dynamo action can occur. 

The asymptotics of helical dynamos in an infinite cylinder have been gener- 
alised to helical dynamos in a sphere (GP). A spherical single- roll dynamo can 
be portrayed as a cylindrical helical flow with its two ends joined, bent into a 
donut, and deformed to fill a spherical 'container'. The theory predicts that at 
large i?ni the dominant mechanism of field generation is of Ponamarenko type 
and, indeed, that the same asymptotic scalings hold. However, full numerical 
studies of the dynamo problem in spheres have mostly tested small or interme- 
diate i?ni, primarily to establish the critical i?ni necessary for dynamo action 
(Dudley & James 1989, Forest et al. 2002, Moss 2008). The very large Rm 
regime has not received the same attention, partly because it remains a chal- 
lenging numerical task. Consequently, the GP asymptotic theory has not yet 
been numerically confirmed. This is the main project of our paper 

We undertake an analysis of single-roll helical dynamos in a sphere, sur- 
rounded by an insulator, at very large magnetic Reynolds numbers. This is 
accomplished numerically by the solution of the magnetic induction equation 
with suitable boundary conditions for two representative single-roll flows. The 
method of solution is presented in Ivers and Phillips (2003, 2008) and consists of 
approximating the problem as a large-scale algebraic eigenvalue equation. The 
growth rates and eigenfunctions so obtained are subsequently compared with 
the predictions of the GP asymptotic theory. We find excellent agreement in 
the regime Rm > 10* for the growth rates and angular frequency, and good 
agreement for the magnetic field structure when R^ > 10^. In particular, the 
magnetic field structure clearly localises upon a specific stream-surface, in con- 
trast to intermediate R^ where the field is somewhat disordered. This confirms 
that spherical single-roll dynamos are indeed of the Ponamerenko type: no other 
growing modes were discovered. This is a point that can deepen our understand- 
ing of spherical dynamos more generally and aid the analysis of more complex 
spherical fiows. In particular, it is an entry point into the study of multiple rolls, 
whose magnetic field generation will also be influenced by the Gailitis dynamo 
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mechanism (Gailitis 1970, 1993, 1995, Moss 2006). 

The layout of the paper is as follows. In Section 2, the formal dynamo prob- 
lem is stated — its governing equations, parameters, and boundary conditions 
— and the two flows we examine arc presented. A brief summary of the GP 
asymptotic results at large i?m follows in Section 3, while their derivation is 
given in Appendix A, in the Supplemental Material, alongside a method to ob- 
tain higher order terms. Our results and a comparison of the two approaches 
are given in Section 4, and conclusions drawn in Section 5. 

2. Governing equations and setup 

2.1. Problem formulation 

Consider a sphere of conducting fluid V with radius a and uniform magnetic 

diffusivity rj surrounded by an insulator V^. Suppose that the fluid is undergoing 
time-steady incompressible motions according to the velocity v. Consequently, 
the magnetic field in the conducting fluid is governed by the non-dimensionalised 
induction equation, 

drB = V^B + R^V X (v X B) (2.1) 

where the magnetic Reynolds number Rm = Va/r] is defined in terms of a typical 
velocity V, the radius a, and rj. The time r is scaled on the magnetic diffusion 
time C^/ri and space by a. The magnetic field B is solenoidal everywhere, 

V-B = 0. (2.2) 

Because the flow is steady, the magnetic field can be expressed as a linear 
superposition of time-separable solutions of the form 

B(r,r)=B(r)e^^ (2.3) 

possibly with polynomial factors of time in degenerate cases. In addition, we 
must supply suitable boundary conditions 

[B]s = , V X B = inV, B^O as r ^ oo, (2.4) 

where S is the surface of the sphere. This leads to an eigenvalue problem for 
the (complex) growth rate A and the associated eigenfunction. For a given fiow 
V, the growth rate A is a function of the magnetic Reynolds number R^. When 
Re A > 0, the fiow acts as a kinematic dynamo, i.e. the non-magnetic state 
B = is unstable to magnetic perturbations. 

2.2. Representation of the helical flow 

We use spherical coordinates whereby the radius, polar angle, and azimuthal 
angle are denoted by r, 0, and (j), with their accompanying unit vectors given 
by If, le, and 1^, respectively. An axisymmetric single-roll helical flows may 
be represented by 

v = aV„, + W{r,e)rsmei^, (2.5) 



3 



^ and -W^ 




Figure 1: A representation of the meridional circulation shared by both flow 1 and 2. The 
coloured isocontours of the meridional stream function ^ are plotted. The negative of this 
function is also the azimuthal angular speed of flow 1, as Wi = — Si/. 



where V„i is the (scaled) meridional velocity, W is the azimuthal angular speed, 
and (7 is a parameter measuring the relative strengths of the meridional and 
azimuthal motion. The meridional flow Vm can be written in terms of a stream 
function ^ by 

Vn. = -^^lr + ^^le=V(/)X V^-^-Vx — (2.6) 
sm tl r sm r sm 

The streamlines of V,„ in a meridional plane are the level contours of 5* and 
circle a local minimum (maximum) of ^ in the clockwise (counter-clockwise) 
direction. We also introduce the 'unsealed' meridional velocity Vm = (tV,„ 
and stream function -0 = ct^, which are more convenient when describing the 
asymptotic theory. 

In this paper, we examine two representative single-roll flows, vi and V2. 
The two flows share the same meridional velocity but differ in the azimuthal 
component Wr s'm6. We set the stream function according to 

*(r, e) = -r sin vrr sin^ 9 , (2.7) 

but set 

Wi=r sin nr sin^ 9, W2 = . (2.8) 
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Figure 2: A representation of the two azimuthal flows considered. The left panel shows the 
azimuthal velocity of flow 1, Wi rsind, while the right panel shows that of flow 2, W2 rsind. 



Flow 1 therefore possesses the restricted form Wi = Wi{^) = — Vf, while flow 
2 has a more general form. The restricted form of the local angular velocity 
simplifies the asymptotic theory substantially. 

In Figs 1 and 2 the two flows are represented. Figure 1 shows the isocontours 
of the meridional stream function common to both flows. As the figure indi- 
cates, fluid circulates in a clockwise direction around a central stagnation point 
located at r « 0.646. In the case of fiow 1, for which Wi — — ^, the isocon- 
tours of the azimuthal rotation coincide with the isocontours of the meridional 
rotation. Therefore a fluid element upon a given streamsurface will not only 
circulate at a constant speed in the meridional plane it will also travel at a 
constant angular speed in the azimuthal direction. It is then easier to flnd flow 
trajectories that form closed loops. 

This is not the case for the more complicated flow 2, whose azimuthal angular 
speed exhibits a purely radial profile. Thus surfaces of constant rotation are 
spherical shells. This means a fluid element in flow 2 will experience different 
azimuthal angular velocities as it traverses a meridional streamline. The actual 
azimuthal velocities, WrsinO, of both flows are plotted in Fig. 2. 

The full numerical method we employ requires that the velocities are split 
into toroidal and poloidal components and expanded in spherical harmonics: 

m,n m,n 

where the toroidal components are given by t™ = V x {t^Y™r) and the poloidal 
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components by s™ = V x V x (s^Y^r), where is a spherical harmonic (see 
Ivers & Phillips 2003, 2008). 

Flow 1 has the poloidal-toroidal spectral decomposition vi = asj + t5 +t3, 
in which the radial functions are 

n sinTrr r, 4 , . r, 2 , . ,^ 
s\ = — , t? = ^ sinTrr , K = ^ sinTrr , (2.9) 

and the spherical harmonics are = a/S cos 9 and Y^ = cos 0(5 cos^ 9—2>). 

The second flow V2 is the single-roll flow of Dudley & James (1989). This 
flow has the poloidal-toroidal spectral form V2 = crsj + tj, where 

n sinTrr n . , 

s\ = —j^, t?= sinTrr. (2.10) 

Although vi has the simpler restricted form of the local angular velocity, its 
spherical harmonic representation is actually more complicated than V2. 



3. The GP asymptotic theory 

This section suinniariscs the main features and results of the GP asymptotic 
theory, undertaken at largo i?^. Here we give, without proof, the leading order 

— 1/2 

terms for the growth rates at 0{Rra ) and the leading order contribution to 
the magnetic eigenfunction. The derivations of these expressions are placed 
in the Appendix A of the Supplemental Material for reference. There we also 
present a method whereby higher order terms in the asymptotic expansion can 
be calculated. 

3.1. Toroidal co-ordinates for axisymmetric helical flows 

The structure of the hehcal flow (2.5), i.e. the topology of its streamlines and 
its differential rotation, can be exploited by a toroidal coordinate system which 
simplifles the advection operator. The first two coordinates are determined by 
the meridional flow. The stream function tp is one coordinate (note that we 
use the unsealed stream function). The second is an angle coordinate i? deflned 
as follows. If T is the period for a fluid particle to traverse once the closed 
streamline tjj = ipo, then T = T{tjj) and 

2Tr , ^di ^'Vm-dr ^rdO ^dr 

di) = —dt = n— = n " =n — = (3.11) 

T q ve Vr 

where f2 = 2'k/T is the angular frequency, d£ is the element of arc- length trav- 
elled in a time dt and q = |vni|. Thus dd/dt = ft is constant on streamlines. 
For the two flows we consider we flx ■!? = on the s-axis {6 = tt/2), since the 
stagnation points of their meridional parts occur there. Clearly i? changes by 
2Tr in one full traversal of the closed streamline. We denote by an overbar, or 
pair of angle brackets, the average around the streamline tp = tpo, defined by 



/ = (/) = 




fd^, 
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where / is a function of the meridional coordinates "0 and and the integration 
is around the streamhne tp = V'o- We can thus take any quantity F{^,t{j) and 
compute a mean component independent of i?, i.e. F = F{%p), and a 'fluctuating' 
component, F{tjj,i')) = F{ijj,'d) — F. 

Having specified these details, we replace the azimuthal angle ^ by a third 
coordinate ( defined by 

C(v,^,<^) = </>-^(V,t?), z = -^ j^w{^p,r)dr . (3.12) 

The level surfaces of ^ are hence distorted azimuthal planes 

The coordinate system () naturally gives rise to the two right-handed 

vector bases, (VV', Vi?, VC) and (f^,f^,f(^), where fj/, = dr/dtp, = dr/d'd, 
fQ = dv / dC, and r is the position vector. It is a useful shorthand to also denote 
the coordinates by with indices i = 1, 2, 3, and the two bases by V^* and fj. 
The bases are reciprocal, fj • VV'-' = , and related by 

= JVi? X VC, U = JVC X , fc = JVip X Vi? , (3.13) 

where the Jacobian J of the transformation to {ip, Q is given by 

J = X ftf • = (V^/- X Vi? • VC)"^ = 0"^ (3.14) 

Using the properties of the flux function, equations (3.13) may be simplified to 

i^=n-'^V'&xV(p + fc:d^Z, = Q-^V(?i X W + fc^^^, fc=rsin^l0. 

(3.15) 

Finally, using the properties of reciprocal bases the velocity and the magnetic 
field may be written as 

V = n(V')ftf + W(^)fc , B = B^i^ + B^U + Bcfc • (3.16) 

The advection operator is expressed as 

Dt = dt + v-W = dt + n{i;)di) + W{ij)dc . (3.17) 

Its sole dependence on ^ is essential in the asymptotic theory. 

The magnetic induction equation scaled on the turn-over timescale t = ri?m 

is 

a^B = V X (v X B) + £^V^B , (3.18) 

where s = K^^^. The contravariant components of (3.18) with respect to the 
new coordinates are 

DtB^ = £*VV' • V^B (3.19) 
DtB^ - n'{ip)B^ = e'^V'd ■ V^B (3.20) 
A^c - W'{iIj)B^ = e^VC • V^B . (3.21) 

The primes indicate derivatives with respect ip. As noted above, regeneration 
of the magnetic field component B^ is solely due to diffusion, but regeneration 
of B^ and B^ is partly due to distortion of B^ by meridional and azimuthal 
differential rotation, respectively. 
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3.2. Asymptotic scalings 

The GP theory is developed in the large -Rm regime, i.e. as £ ^ 0. The lead- 
ing order solution is decomposed into dynamo modes of the form B^, B^,B(^ cx 
^im(+iM+{p+iuj)t ^ Thus p+ito = X/Rm- The constants m and k must be integers 
for solutions single-valued in ( and t?. The following scalings are subsequently 
adopted (Ruzmaikin et al. 1988) 

e^B^ r^B^r^BQ, p = 0{e^) , m,k = 0(1) , (3.22) 

moreover, it is assumed that modes localise upon a stream surface ^/i = t/jq in a 
layer of thickness 0{e). Thus d- and ^-derivatives are 0(1), but i/'-derivatives 
are 0{e~^). This suggests a new variable T defined through 

xP = xPo + eT. (3.23) 

so that T-derivatives are 0(1). 

The magnetic solution and growth rates are subsequently expanded in pow- 
ers of e, and a particular streamline V'o is chosen and equilibrium quantities 
depending on ^ expanded in Taylor series about this streamline. These are 
substituted into the governing equations and we collect terms order by order 
solving each set of equations as we go. These details arc summarised in Ap- 
pendix A of the Supplemental Material. In brief, solvability of the zeroth order 
equations forces a given dynamo mode to localise on its critical stream surface, 
where it will be convected with the flow. Solvability of the first order equations, 
yields a 'resonance condition' which selects which stream surface (or surfaces) 
can actually harbour such dynamo modes. Solvability of the second order equa- 
tions provides the leading order T-structure of the modes and the leading order 
growth rate. 



3.3. The leading order asymptotic solution 

In the following expressions we transform from ^ to \E' and thus make explicit 
the dependence on the tuning parameter a. 

3.3.1. Eigenf unctions 

To dominant order spatially, dynamo modes take the form 

B (X a£»„(r/K)e*'=''+™'^+*"'*+f* (3.24) 

where Dn{z) is the parabolic cylinder function of order n (Abramowitz & Stegun 
1972), with 

S. = an'j^ + W'j^ = l^\„, + ^(^W + W')j rsin^l^j, (3.25) 



and 
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where n = akfl + mW, 79 = (|V^'oP), and a prime now indicates differentia- 
tion with respect to ^1/. The parabohc cyhnder functions impart a Gaussian- hke 
structure about the resonant curve, with spatial oscillations of rapidly diminish- 
ing amplitude as distance T from the curve increases. Higher n modes display 
more complex spatially varying behaviour within the envelope of the stream 
surface localization. 

3.3.2. Growth rates 

The real part of the growth rate is 

p = ±eWmb + mJlM - (n + l)e'^\kan- + mCllo + 0{e'') (3.26) 
and the angular frequencies are 

CO = -akCto - mWo ± vl^^T+^^^^clP^ sgn [(fc/Zf, -|- rh4i^)Q!g] 

- (n + i)£2y^|am^' + mW"|7o sgn (akn^ + mW'i^ + 0{e^) , (3.27) 

where the quantities appearing are evaluated upon the resonant stream surface 
\1/ = "^o- The n = mode is the fastest growing magnetic field mode to dominant 
order for any m and k. The geometric fj. terms are 

= o-V* • (Vi? • Vfi») , = V-S- • (VC • Vi?). 

In addition GP add terms of higher order e'^ to the expression for p arguing 
that these are the most important at this order. The additional terms are 

The /3's, which are independent of a, are given by 

= I Vi?|2 - ^ (V* • S/df (3.28) 
7o 

/3m = |VCP-^(V*-VC)' (3.29) 
7o 

Pmk = Vi? • VC - ^ (V* • Vi?) (V* • VC) , (3.30) 

7o 

3.3.3. Resonance condition 

The expression for the mode frequency Eq. (3.27), shows that dynamo modes 
are (to leading order) 'advected' by the flow at the streamline upon which they 
localise. That is to say, w « — B • Vv. It follows that the resonant streamsurface 
is a magnetic critical layer, and modes in its immeditae vicinity will naturally 
exhibit rapid spatial oscillations, though these are regularised by the (small) 
magnetic diffusion upon the critical streamsurface itself. At large these 
oscillations are crucial to dynamo action, because they provide sufficiently steep 
spatial gradients for the small resistivity to work efficiently, and regenerate -B^. 
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They hence close the dynamo loop begun by the differential rotation across the 
layer. 

We have not yet stated which streamline a given dynamo mode will prefer, 
upon which the resonance condition holds. This condition may be written as 

fco-n; + mW'„ = , o-^n;' = kafl'^ + mW" ^ . (3.31) 

In general, we find that 11" < 0, which indicates that the resonant streamline 
corresponds to the maximal helical gradient of the magnetic mode. But Equa- 
tion (3.31)(a) is also the condition for the closure of the magnetic field lines on 
the surface ip = ipo to leading order, as the following argument shows. 

Since dr = i^dtp + f^d-d + f^rfC, the equation for the magnetic field lines, 
B X dr = 0, reduces to 

dtp _ d-d _ dC 
B^ 

For the magnetic field Eqs (3.24)-(3.25), the field lines to leading order are 

w 

where (\E'o, i^oj Co) is a given point on the field line. The magnetic field line is 
closed if there are integers k, m such that d — do = —2iTk and Q — Co = 27rm, 
which give the resonance condition (3.31)(a). So, unsurprisingly, a resonant 
surface also corresponds to the spatial localisation for which a magnetic mode 
reinforces itself. The resonance condition also ensures the e°-magnetic field is 
solenoidal. 

3.4^. Discussion 

In practice, it is simplest to stipulate the ratio (m/k) and the streamline 
^ and then compute the a necessary for this ^ to be resonant from (3.31). 
Therefore, 

^ = (^) =^(Wfc,*o), (3.32) 

and a may be interpreted as an adjustable 'tuning' parameter, permitting modes 
on any streamline we choose. 

However, it is also instructive to examine how the dynamo modes, and their 
resonant streamlines change, as a varies. The parameter a controls the geometry 
of the helical flow by establishing the ratio of the meridional circulation's speed 
against the azimuthal rotation, and hence directly influences the dynamo action. 

Consider flow 2. According to Eq. (3.32) each choice of a and m/k will give 
a single resonant curve Now flx the pitch m/k of the dynamo modes under 
consideration; this means that as we vary a we also vary the resonant streamsur- 
face. On the other hand, for flow 2, it can be shown that, as ^ varies between 
and its minimum value, the quantity Wg/^'o varies monotonically between two 
nearby constants, the smaller associated with the outermost streamline upon the 
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spherical boundary, and the larger with the stagnation point. From (3.32), it 
then follows that there exists only a (narrow) interval of cr for which a resonance 
on any streamline is possible. Each choice of m/fc furnishes a different interval 
of a but none of these overlap. A subintcrval of each may permit magnetic 
growth. Therefore as we vary a, and consequently modify the flow geometry, 
we encounter discrete 'windows', or 'resonance intervals', of magnetic field gen- 
eration. Note that dynamo action is not possible in every interval, in particular 
for very small and very large cr. These require large m or k which violates the 
scaling assumptions of the asymptotic theory, Eq. 3.22. In any case. Proctor's 
modification of the toroidal anti-dynamo theorem (Proctor 2004) suggests that 
there can be no magnetic growth for a < e^^^, i.e. for flows almost entirely 
azimuthal. In the other limit, a large, which corresponds to a flow dominated 
by the meridional component, things arc less clear. For certain stream func- 
tions dynamo action appears possible in the complete absence of the azimuthal 
motion (cr — )• oo), though the relationship between the flow and the boundary 
is crucial (Moss 2008). 

In contrast, the simpler flow 1 admits dynamo action for a very wide range of 
cr; moreover, multiple m/k modes may grow concurrently. In other words, the cr 
resonance intervals overlap siibstantially. Plainly, a simpler flow, in which both 
the azimuthal and meridional motion are closely related, is the more propitious 
for magnetic generation. This follows from the fact that the fluid trajectories 
can more easily join the magnetic field lines they convcct into closed loops. In 
circumstances where the profiles of \l/ and W are dissimilar (flow 2), this is 
more diflicult to do, and can only occur when their relative magnitudes are 
tuned appropriately (by cr). 

Considering how vital it is to for the flow to close magnetic field lines in the 
large Rm limit, it is natural to enquire into the (possibly deleterious) influence 
of small velocity fluctuations superimposed upon the mean helical motion. Re- 
cent work on the cylindrical Ponomarenko dynamo shows that magnetic growth 
persists when the amplitudes of the helical flow has a small time-dependent 
(fluctuating) part. Dynamo action even can occur when the meridional and az- 
imuthal fluctuations are slightly different functions of time, forcing the resonant 
curve to also change with time (Peyrot et al. 2007, 2008). Similar behaviour 
undoubtedly carries over to the spherical single roll dynamos we consider. How- 
ever, when the small velocity fluctuation is not only a function of time, but of 
space as well, dynamo action will most likely suffer. In such a flow the fluid 
trajectories will not normally close. Thus, like flow 2, the magnetic field lines 
they transport will not normally close, and the resonance condition will be more 
difficult to satisfy. Dynamo activity may still be possible in the limit of small 
fiuctuation amplitude, as then the fluid trajectories may not deviate beyond 
the magnetic localisation, but this probably can only be checked with numerical 
simulations. 
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4. Results 

We present results for each of the flows Vi and V2, corresponding to a repre- 
sentative configuration of the parameters for the same resonant stream surface. 
For both vi and V2 this resonant curve is \l/o = 7r/2) « —0.20287 with 

Tg = 0.93. The resonance is ensured for given m and k by setting the tuning 
parameter a according to (3.32). The resonance conditions for Vi and V2 may 
be expressed as cr = ai{m, k, and a = 0-2 (m, k, 5") respectively. We find that 

cri(l,l,*o) « 0.1373, fTi(2,l,*o) - 0.2747, 0-2(1, 1, *o) ~ 0.2050, (4.33) 

to 4 significant figures. Moreover, for the flows we examine there is no degener- 
acy in fT, i.e. for a given a there is only one possible set of (m, k, \l/), and hence 
only one resonant curve for a given flow. 

The different times t, t of the asymptotic and numerical results are related 
by t/r = i?ni- Thus to compare the asymptotic and numerical results, the 
numerical growth rates are divided by Rm, i.e. p = ReX/Rm, w = ImA/iim- 
Moreover, asymptotic and numerical modes must be correctly matched. There 
is no difficulty with the azimuthal wavenumber to, since it coincides in the 
asymptotic and numerical results for vi, and also for V2, when ( is decomposed 
and the factor e*™"^ is extracted from the eigenfunction. We assumed that, 
once the resonant curve and the wavcnumbcrs to, k are chosen, which sets 
the tuning parameter a{m,k,^), the collection of modes determined by the 
numerical eigcnproblcm correspond to the various asymptotic n modes. The 
strongest growing exact (numerical) mode was identified with n = 0. 

4.1. Numerical methods 

The calculation of the asymptotic solution for a given flow and its rendering 
in spherical coordinates is not straightforward. This called for various analytical 
tricks and numerical techniques, a full explanation of which wc give in Appendix 
B in the Supplemental Material. Most of the effort lay in computing the and 
(3 coefficients and also in determining the ^ derivatives of SI and W at ^ = ^o- 
These quantities, in fact, can all be expressed as contour integrals of various 
kinds over all or part of the streamline \l/o. When these integrals are closed they 
can be numerically approximated with excellent accuracy. The chief numerical 
parameter here is the number of (r, 9) points used to discretise the resonant 
streamline 4*0; this was denoted by K (see Appendix B in the Supplemental 
Material). 

On the other hand, the full dynamo problem presents a considerable numer- 
ical challenge, especially in the large i?m regime when the magnetic structure 
exhibits small-scale localised variation. In this limit extreme resolution is re- 
quired to properly capture the dynamo modes, which leads to the eigensolution 
of enormous, albeit banded, matrices. A full explanation of the techniques em- 
ployed to approximate the magnetic induction operator with such matrices can 
be found in Ivers & Phillips (2003, 2008). We summarise the approach below. 
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Re A 


ImA 


J\N 


20 


30 


40 


20 


30 


40 


200 
400 
800 


687.7 
688.2 
688.2 


688.2 
687.5 
687.5 


688.3 
687.5 
687.5 


16390.6 
16390.6 
16390.6 


16391.6 
16391.9 
16392.0 


16392.8 
16392.1 
16392.3 



Table 1: Convergence of the dominant mode's eigenvalue A with the numerical parameters J 
and N. The number of radial grid points is J + 1, the number of spherical harmonic functions 
is 2{N — m + 1). The growth rates are for V2 at Rm = 10® and m = 1. 

Our numerical method uses a hybrid version of the spectral equations of 
James (1974), which are in a similar form to the poloidal-toroidal spectral equa- 
tions derived by BuUard & Gellman (1954). The magnetic field and velocity field 
are expanded in vector spherical harmonics when they appear in the advection 
term, but otherwise are decomposed into toroidal-poloidal components and then 
expanded in scalar spherical harmonics. Doing so reqiiires us to compute fewer 
coupling integrals. The spectral expansion is truncated at some large order N, 
so that for a given azimuthal mode number m there are 2(A'' — m + 1) spherical 
harmonic functions. 

The radial dependence is discretised using fourth-order finite-differences over 
a uniform grid. The number of radial points is denoted by J -I- 1 . A centred- 
finite difference formula were used at interior points and one-sided formulas at 
the boundaries. 

Truncation in the number of harmonics and in radius converts the problem 
into a set of linear equations for (2 J -|- 1)(A^ — m -|- 1) coefficients plus the 
growth rate A. We solve this algebraic eigenvalue problem by inverse iteration 
and the implictly restarted Arnoldi method using ARPACK (Sorenscn 1992). 
The Arnoldi method is particularly helpful in identifying the mode of fastest 
growth at large i?m because the eigenvalues in this limit bunch together in the 
complex plane: the ratio of the real part of the growth rate to the imaginary 
part is 0{Rni ) as indicated by the asymptotic theory (see also Table 2). 
Consequently, inverse iteration has difficulty in converging to the eigenvalue of 
largest real part without a good estimate of the true eigenvalue. 

At very large magnetic Reynolds number {Rm ^ 10^)? a- converged solution 
(with respect to resolution) requires extremely large J and N. The largest 
matrix computed was 64,040 x 64,040 for J = 800, TV = 40, and m = 1. 
Convergence of the eigenvalue A of largest real part with respect to J and N 
is shown at i?m ~ 10^ for V2 in Tabic 1. The higher order modes require even 
greater truncation levels, as they exhibit steeper spatial gradients (being more 
strongly localised). As Ry^ is increased further (i^m > 5x 10^) difficulties are 
encountered because the truncation needed (and consequently the size of the 
matrices generated) become prohibitive. 

4-2. Growth rates 

We present below the growth rates of the n = 0,1,2 magnetic modes for 
k = 1 and m = 1 in the two spherical helical dynamos we considered. 
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Re Ao 


Im Ao 


ReAi 


Im Ai 


ReA2 


Im A2 


500 


17.8 


18.2 


-60.3 


17.8 






1,000 


38.0 


74.7 


-48.0 


66.1 


-642.4 


113.6 


2,000 


68.9 


206.3 


-25.0 


187.4 


-384.8 


196.7 


3,000 


93.6 


347.2 


-5.2 


320.8 


-375.0 


322.2 


5,000 


133.2 


641.3 


27.4 


602.2 


-348.1 


592.7 


10,000 


203.8 


1408.9 


87.8 


1344.3 


-282.1 


1303.8 


20,000 


293.9 


3006.6 


178.5 


2878.7 


-211.5 


2747.8 


30,000 


363.0 


4645.6 


224.4 


4423.7 


-237.0 


4339.7 


50,000 


477.1 


7969.4 


250.7 


7697.2 


-86.1 


7546.4 


100,000 


687.5 


16392.3 


348.9 


16000.0 


-19.2 


15566.2 


200,000 


982.4 


33432.1 


486.2 


32867.1 


-57.2 


32299.4 



Table 2: The growth rates of the leading modes as computed by the numerical eigenproblem 
for V2, m = 1 at different Rm- The subscript on A indicates the mode number n. Numerical 
truncation levels are AT = 40 and J = 800 

Table 2 shows the growth rates of the n = 0, 1, 2 modes, as computed by the 
numerical eigenproblem for V2 at the truncation levels = 40 and J = 800. 
For sufficiently large Rm, there exist two growing dynamo modes corresponding 

1 /2 

to n = and n — 1. Note the scaling ImA/RcA ^ Rni predicted by the 
asymptotic theory for large i?m (see Equation (3.22)). Moreover, the leading 
modes possess growth rates whose imaginary parts asymptote to a common 
value (see Equation (3.27)). This characteristic clustering of the eigenvalues 
in spherical helical dynamos explains the difficulty that algebraic eigensolvers 
encounter when separating the eigenvalues at large Rm- In this regime a partial 
eigensolver, such as the implicitly restarted Arnoldi method, is invaluable (see 
Latter and Ivers 2004). 

In Figs 3 and 4 we directly compare the predictions of the asymptotic theory 
and the numerical computations of the full eigenproblem for vi. Here is pre- 
sented the numerical growth rate Re A, and associated angular frequency Im A, 
of the leading modes as a function of Rm. Alongside these data points, we 
plot the asymptotic growth rates pi and frequencies LUi , calculated with (a) the 
higher order terms of GP included (the solid lines), (b) the asymptotic theory 
correct to order (the dashed line). The truncation levels for the asymptotic 
values arc K = 400. Note that both growth rates are scaled on the turnover 
time. This makes clear that the numerical growth rate goes to zero as Rm -> 00. 
Helical dynamos are 'slow', as expected (Childress and Gilbert 1995, GP). Of 
the three modes shown, only the n — mode is a dynamo, which is active above 
a critical magnetic Reynolds number of i?m ~ 416 and achieves its maximum 
positive growth rate at i?m ~ 1500. In Fig. 5 we plot the growth rates and 
frequencies of the leading two modes of the V2 dynamo, for m = k = 1. Both 
of these modes may grow for sufficiently large i?m- 

Plainly, there is excellent agreement between the n = asymptotic growth 
rates and the numerical growth rates when > 10** for both vi and V2. The an- 
gular frequencies agree at smaller i?ni- Fig. 3 also shows that the additional terms 
of GP at higher order improve the accuracy of the asymptotics markedly, in 
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Figure 3: The growth rates and angular frequencies of the two leading modes of the vi flow. 
These are characterised by m = 1, fc = 1, and n = or n = 1. The solid line represents 
the asymptotic expressions with the extra terms of Gilbert & Ponty (2000), the dashed line 
represents the asymptotic expression correct to order e^, the points on the dotted line represent 
the full numerical eigensolution. 




Figure 4: The growth rate p2 and angular frequency 1^2 for the m = I, k = 1, n = 2 mode of 
VI . 
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Figure 5: The growth rates and angular frequencies of the leading two modes of the V2 
dynamo. Here m = 1, fc = 1 and n = or n = 1. 

comparison with the theory up to (the dashed line). This agreement strongly 
supports the asymptotic theory. It also indicates that the identification of the 
numerical modes with the asymptotic modes is correct. Finally, this shows that 
these flows only admit Ponomarenko-type dynamos in the large R^n regime. 

4-.3. Magnetic field structure 

Figures 6-8 show the real parts of the magnetic eigenfunctions for different 
values of R^n and n upon the Vi flow. These choices reveal the salient physical 
and asymptotic features of these modes. The morphologies of the V2 dynamos 
are much the same and are omitted in the interests of space. 

In Fig. 6 we plot the magnetic field components, B^, Bg, i?^, of the m = 1, 
fc = 1, n = mode. Here cr = cri{l, l,^o) and we set the magnetic Reynolds 
number to i?m = 10^. The top three panels present the numerical eigenso- 
lutions, while the bottom three panels show the asymptotic approximations. 
Superimposed upon the first set of figures is the resonant stream curve, which 
helps highlight the localisation of the magnetic field. 

It is apparent from the figures that the asymptotic and numerical magnetic 
fields agree in their dominant features: the position, orientation and shape of 
the local maxima and minima. The localisation of the field to the resonant 
streamline is readily observed, with marked fiux expulsion inside and outside 
the resonant streamline, as expected. In addition, the k = 1 nature of the 
field is clear from its variation around the streamline, especially in B^. Note, 
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Figure 6: Plots of the magnetic eigenfunctions: the real parts of Br, Be, B^, for the m = 1, 
fc = 1, n = mode upon the vi flow. Here a = o"i(l, 1, 'to) and Rm = 10^. The three upper 
panels present the numerical eigensolution, while the lower three panels present the asymptotic 
eigensolution. Both numerical and asymptotic eigenfunctions have been normalised so that 
jB| = 1, and each component has been plotted with 10 contours of equal increments. 
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however, that there is a shght offset outwards away from the resonant strcamhne 
in the upper numerical eigenfunctions. This offset in combination with the 
steep spatial gradients means the relative error between the asymptotic and 
numerical eigenfunctions is larger than expected. The relative error, defined by 
||B„um _ B'^''y™||2/||B™™||2, is 0.334 for this = 10^ case. 

In Fig. 7 we show the same mode but for larger magnetic Reynolds number: 
i?,n = 5 X 10"''. The resonant stream curve has been omitted. The magnetic 
field is more localised and intense, with the offset of the maxima and minima 
from the resonant streamline substantially reduced. Now the magnitudes of the 
numerical and asymptotic eigenhmction, in addition to the orientation of the 
magnetic features, are in good agreement. The relative error is reduced to 0.285. 

The more involved n = 1 eigenmode is plotted in Figure 8. The other 
parameters are kept the same and = 5 x 10^. This mode decays slowly 
with time as shown in Figure 3. Its structure is more complicated due to the 
variation with ip under the gaussian envelope of £>i(T//t). In particular, the 
field components vanish on the resonant streamline. The agreement between 
the asymptotic and numerical eigenfunctions is good but not as striking as in 
the n = case. In particular, the relative magnitudes of the maxima and 
minima are not in agreement, though the general shape of the structures are. 
This is possibly due to the fact that the mode has not fully converged to its 
asymptotic form. On the other hand, its greater spatial variation may be taxing 
the resolution of our numerical scheme, leading to errors. 

5. Conclusions 

In this paper we have compared the asymptotic theory of Gilbert & Ponty 
(2000) for axisymmetric roll dynamos in a sphere to the numerically computed 
results of the exact dynamo problem for two simple flows, with azimuthal com- 
ponents of the special form = rsmdW{ip) and of a more general form. In 
the regime Rm ^ 10^ excellent agreement is obtained between the asymptotic 

— 1/2 

theory to 0(i?m ) and the numerical results for the growth rate and angular 
frequency. For the magnetic field the agreement between the asymptotic theory 
at leading order and the numerical results is good if i?m — 10^ and excellent 
if Rm = 5x 10^. The asymptotic formulas for the growth rate and the angular 
frequency have been extended in Appendix A in the Supplemental Material. 

Only the simplest class of axisymmetric roll dynamos have been considered: 
those which consist of a single-roll flow with a single resonant streamline. The 
magnetic fleld in these dynamos is localised to the resonant stream surface and 
can interact only with itself. Further work is required on more complicated 
spherical roll flows: such as those those with a single-roll but more than one 
resonant streamline, and those with several rolls. These flows offer the possi- 
bility of interaction between magnetic fields localised to separate regions of the 
fiow. This may produce interacting modes of non-Ponomarenko type, e.g. Gaili- 
tis type modes, alongside the Ponomarenko type modes (Gailitis 1970, 1993, 
1995, Moss 2006). A related question, which arises from the localised nature 
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of the Ponomarcnko modes, is whether they depend on the magnetic boundary 
conditions at the surface of the conducting fluid. 

Of key interest to both laboratory and astrophysical apphcations are (a) the 
nonhnear saturation of such modes and (b) their relationship to a backgroimd 
of small-scale velocity fluctuations. How will the dynamo modes back react 
on the helical flow which generated them? Will these modes build up magnetic 
torques which stifle the meridional motion, or will a more complicated dynamical 
interplay arise? If fluctuations are indeed present, for what size amplitudes, 
and for what correlation times and lengths, will they succesfuUy impede the 
satisfaction of the resonance condition which is so crucial for the formation of 
Ponomarenko modes? Recent work in cylindrical geometry for time-dependent 
fluctuations shows that Ponamarenko-type dynamo action can, in fact, survive 
in certain cases (Peyrot et al. 2007, 2008). But this need not be the case at all 
for flows exhibiting small-scale variations in space, and it is the latter situation 
which is probably most relevant in applications. 
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Spherical single-roll dynamos at 
large magnetic Reynolds number: 
Supplemental material 

Appendix A. Derivation of cisymptotic expressions at leirge Rm 

In this appendix we briefly derive tlie leading order results of the GP asymp- 
totic theory presented in Section 3 of the main manuscript 'Spherical single- roll 
dynamocs at large magnetic Reynolds number'. Once these are given we sketch 
out a technique whereby higher order terms may be calculated and give the solu- 
tion to order . Note that reference labels for equations in the main manuscript 
are not preceded by either an 'A' or 'B'. 

Appendix A.l. Preliminaries 

Appendix A. 1.1. The magnetic diffusion term in toroidal coordinates 

We take as our starting point the magnetic induction equation in the toroidal 
coordinate system {ijj, d, C), Eqs (3.19)-(3.21) in the main manuscript. In order 
to progress, the diffusion terms on the right sides need to be decomposed into 
their component parts. 

Using the summation convention the magnetic field can be written as B = 
B*£j, the gradient operator as V = (VV'')^^, where di — d^i, and the diffusion 
term in the magnetic induction equation as V^B = (V^5-')fj + 2V-B-' • Vfj + 
B^V^ij. Thus, since VB^ = {dkB^)V'tjj'' , the covariant components of this term 
are 

Vtjj' ■ V^B = (V^B*) + 2(afeB^)VV''= • Vfj ■ VV'* + B^V^)' ■ V^f,- . 

Ibur of the 27 terms Vtp'' ■ Vfj ■ Vtp^ and three of the 9 terms Vii)' ■ V^fy 
vanish identically, since Vf^ = l^l^ — l<^ls means that a - (Vf^) • b = for any 
meridional vectors a, b. We then have 

Vip ■ {Vtp ■ Vfc) = Vip ■ ■ Vfc) = Vi? • (VV' • Vfc) = Vi? • (Vi? • Vfc) = . 

Also V^f^; = implies VV^ • V^f^; = Vi)* • V^f^ = VC • V^f^ = 0. 

Keeping only terms which appear later in the asymptotic analysis and sup- 
pressing the others with dots, the relevant diffusion terms are 

• V^B = (V^ + 2ii,d^ + 2iijds + 2iJkdc + ni)B^ 

+ (2/ia9^ + 2ij,bd^ + 2/Uc9^ + l^d)B^ + 2ngdQBc^ 

(A.l) 

Vd ■ V^B = (V^ -f 2Xad^ + 2\bdi) + 2\^d^ + \d)B^ + 2\gdcBc^ + ... (A.2) 
VC • V^B = 2p^d^B^ + (V^ + 2pbd^)B^ + ..., (A.3) 
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where the coefficients are defined by 

Ma = • (V^ • Vftf ) M5 = • (Vi? • Vf^) Me = V^i- • (VC • Vftf) Atd = V^-V% 

Mg = • (VC • Vfc) = • (VV • Vfv,) = • (Vi? • Vfv,) Mfe = • (VC • Vfv,) 

^, = VV'-V2f^ Aa = Vt? • (VV- • Vftf) Ab = Vt? • (Vi9 • Vfi,) Ac = Vt5 • (VC • Vf^j) 

Ad = Vt?-V% Ag = • (VC • Vfc) p„ = VC • W • Vftf ) p6 = VC • W • Vfc) . 

Apart from /i^, f^h A^ these are ChristofFel symbols. 

The scalar Laplacian is = V • (V^^^i = (VV'* • '^^^)didj + (VV*)5i, 
which, in full, can be expressed as 

= {Vi^ fdl + 2(VVj • Vd)d^di, + (V^ fdl + 2{Vd ■ VC)a^9c 

+ (VC)'52 + 2(VC • Vij)dcd^ + {V^iP)d^ + {W^^)di) + (V'C)5c . 

Furthermore, a number of the geometric coefficients average to zero, 

which can be established using standard vector identities, the divergence theo- 
rem and Stokes' theorem (see GP). 

Appendix A. 1.2. Explicit asymptotic expansions 

We now employ the asymptotic scalings presented in Section 3.2 and the 
order one variable T, and then choose a specific streamsurface tp = tpo around 
which we expand. 

The V-derivatives and gradient operators in Eqs (3.19)-(3.21) become 

d^ = e-^dr, V = £-VV'^ + Vi?a^, + VC^c (A.5) 
and the magnetic field components take the functional forms, 

(A.6) 

We now expand w and p in powers of e with the ordering (3.22), 
w = Wo + + e^W2 + e^ws + e^ws + . . . , p = e^p2 + e^ps + e^p4 + . . . , (A.7) 
and expand W{ijj) and ri(V') in Taylor series about the streamline ip = tpo, 

n{iPo + eY) ^^0 + n'^er + ^n'yr^ + ... (a.s) 

Wi^Po + sT) =Wo + w'^eT + ^W'yr^ + ... , (A.9) 

in which fio = U{ipo), ^'o = f^'(V'o), etc. Assuming the functional dependencies 
of (A.6) and substituting the expansions (A.7)-(A.9) into the advection operator 
(3.17) gives 

Dt = do+ edi + £^d2 + e^d^ + e'^di + ... , (A. 10) 
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where 

d„=p„ + zw„ + — (l^i")a^+imWrj, po=Pi=0. (A.ll) 

We also expand (VV')^ in the diffusion term, 

(W)^ = 70 + £3^71 + £^3^'72 + e'3^^73 + £^^^74 + . . . (A. 12) 

as well as the individual 'diffusion coefficients', for example: 

= /"6,0 + £^M6,1 + 

Finally we expand the magnetic field components, 

61/" ^ b-'l'O + + ■■■■, ^i) = b-m + ehi + • ■ ■ , &c = ^Co + £^Ci + ■ • ■ • (A.13) 
We are now ready to derive the asymptotic equations at the various orders. 

Appendix A. 2. The £° equations 

In this and the following two subsections wc describe the asymptotics to 
order e^. We substitute expansions (A.7)-(A.10), (A. 12) and (A.13) into the 
component equations (3.19)-(3.21), divide (3.19) by and collect terms of like 
order. 

The £°-equations are 

dob.^0 = , dob^ = , do^'co = , (A.14) 

which have the solution, 

6^.0 = ^^0(^)6'*='^, bM = FM{T)e"'\ b^o = F^o{T)e'^\ (A.15) 

where the functions F^q, F^o, -f^o ^■re determined at order and must vanish 
as \T\ — > 00. The constant k is an integer since B is singlc-vahied. Solvability 
of (A. 14) fixes the angular frequency w to leading order for given m and fc, 

wo = -no, (A.16) 

where wc have introduced the advcction frequency function !!(?/') = fcfi + mW . 
The operator do becomes ^o{di) — ik), and hence annihilates any term with the 
i^-dependence e''^''. 

Appendix A. 3. The equations 
The £^-equations are 

dob^i + dib^o = 'i^J.a.odrb-do , dob^ + dibi)o ^ , do^^ci + '^i^CO = . (A.17) 

Their solvability requires 

uji=0, n; = fcfi; + mW'„ = 0. (A.18) 
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The last condition fixes the resonant streamhne tp = tpo, upon which the mag- 
netic field is localised for given m and k. At this streamline the function Il{tp) 
possesses a critical point, and a maximum if 11" < 0, which is the case for the 
simple roll flows we examine. The larger gradients in B{) and i?^ on this surface 
encourage diffusion of these fields and hence replenishment of B^. The oper- 
ator rfi becomes di = Tfl'^{d^ — ik) and hence also annihilates any term with 
■!?-depcndcnce e*'^''. 

The last two equations in (A. 17) can be solved similarly to (A. 14). The first 
equation reduces to, 



\k-9 



which is solvable, since ji^ = Q and the right side then possesses no term with 
the f^-dependence e^'^^ . Thus the magnetic field at order is 

6^1 = F^i(T)e^'='' + G^i(r,t?)e^'='» , = i^^i(T)e^'='» , h^^ = F^^{T)e''^ , 

(A.19) 

where the functions F-^i, F^i, F^i are determined at order and the particular 
integral for equation (A.17)(a) is 

op' 

G^i = -^l2u,o, G^i=0. (A.20) 

Here we have introduced the hat operator ^ defined by 

d^f = f-J, /=0, 

which implies 



/= / {f-f)d^- [ {f-f)d^. 
Jo Jo 

In addition, the properties 

T9 = -l9, 77 = (A.21) 

are easily established. 

Appendix A. 4- The equations 
The e^-equations are 

doh^i + d\h^\ + {di - 7o9r)^V'0 = 2na,odrbn + ^/j-'a^oTdrbM 

+ {2ij.b,od^ + 2imfic,o + IJ'd,o)bi)o + ^imng^obco 

(A.22) 

dnb^2 + dib^i + (^2 - Jodx)bm = f^o^i/.o (A.23) 
dob(;2 + dibf^i + {d2 - 7o5j)6co = W^b^o ■ (A.24) 
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GP included subdominant terms from the Laplacian at this order arguing that 
these are comparable when employing the scalings of Gilbert (1988). However, 
they neglect to include the coordinate Laplacians, V^V) V^i?, V^C which should 
be of the same order. In the present analysis all these terms appear at the correct 
(higher) orders. 

Equations (A.22)-(A.24) are solvable for the field components 6^2, b^2 and 
6^2, if the ?9-dependence of the other terms is not e*'^''. This is true for the terms, 
dib^i, dib^i and since the operator di annihilates e^'^^ . It is also satisfied 

by the terms, 2iJia,odrb^\, '2Hd,obm and 2im^,g^obc,o-, 

since the coefficients iJia,oi 

Hd.o, Mg,o average to zero by (A. 4). The sum of the remaining terms must 
average to zero after multiplication by e~^^'^ . Gonsequently we may write the 
solvability condition for (A.22)-(A.24) as 



_ (A.25) 
where 

S = 7oS^ - |«n"T^ - P2 - , a = fc/if, + m/itc , H" = fcO" + mW 
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(A.26) 

with 7o =< iVV'lo >• These equations determine the functions i^^o, F^q, F(^o 
and hence the magnetic field to leading order. The solutions are of the form 

i^^o = 2/„(r)a , y„(r) = D„ {Y/k) , k = {%/2my/\ (A.27) 

Here a is a constant vector to be determined and D„(2) is the parabolic cylinder 
function of degree n, 

D„(^) = 2-"/2e--V4H„(^/V2), n>0, 

where H„(z) is the Hermite polynomial of degree n (Abramowitz & Stegun 
1972). In order for this solution to satisfy the boundary conditions, i.e. F^o, F^o, -^co 
as |T| — >■ 00, we choose with positive real part. 



= v'|n;'|/7o (i + isgnn;'), 

noting that 7q > 0. The two eigenfunctions which arise from the ambiguous 
sign of K differ only if n is odd and then only in sign. Thus 

The Un are eigenfunctions of S with eigenvalue 



2y„ = ^nVn , ^n = -{n+ l)V\^oho (1 + « sgnH^') - (p2 + iu>2) . 
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Substitution of the ansatz (A. 27) into (A. 25) gives 

/ C„ 2iao \ 
L„a = , L„ = hX, ?n , (A.28) 

u 

which has non-trivial solutions if det L„ = 0. This yields ^„ = or = 2zaof^o: 
i.e. = ±(1 + «sgnao^^o)'\/l^of^ol) ^^'^ determines p2 and W2- We shall ignore 
the solution for ^„ = 0, since its growth rate has a negative real part. The other 



two solutions give for the nth mode 

P2 = TV\^\-{n+l)Vm% (A.29) 

0^2 = tV\^\ sgn(aon'J - (n + ^)V\U^o ^E^K (A.30) 

a=K„,f2;,WlF. (A.31) 



The real and imaginary parts of the growth rate in (A.29) and (A.30) agree 
with Gilbert & Ponty (2000) to order e^. The vector a is determined up to a 
constant factor. 

Equations (A.22)-(A.24) have solutions of the form, 

G^2 = (A.32) 
G^2 = (A.33) 
Gf2 = . (A.34) 

The particular integrals Gy,2, G^2 and G^2 can be determined at this order by 
subtracting from (A.22)-(A.24) their projections on e'*''' and integrating with 
respect to i}. The results are 

2f2' 

OoG^2 = --^tia,oTF^o+%F;i,o+2fla,oF^i+2li'^^^rF^o+{2iao+Jid,o)FM+2im^^ 

(A.35) 

OoG*2 = 7oF^o ■ (A.36) 
i^oGc2 = JoF^o ■ (A.37) 

Appendix A. 5. The asymptotic theory to 0{e^) and beyond 

In this subsection we compute the; solution of the equations and thus 
determine the magnetic field structure to £^ and that Ps = Ws = 0. In so doing 
we sketch out a general method to compute the solution at higher orders. 
The e^-equations are 

+ dib^2 + {d2 - lodr)b-4,i + {ds - 7ira^)6^o = 2(xo,i + Mi,o)9r6^o 
+ 2na,odrb^2 + ^n'^^^Tdrbsi + Ha^^T'^drbm + {2lJ.b,od.& + 2imiXc,o + IJ-d,o)bm 
+ (2/i'(,_„5^ + 2im/i;„ + n'd,o)'^bM + 2imtig,obci + 2im^i'g^Jb^Q (A.38) 



b^2 = F^2{r)e'''^ + G^2(T,{})e'' 
b^2 = Fs2{T)e''''^ + GMT,^)e'' 
b^2 = F^2{r)e^^^ + G^2{T, t?)e''= 
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dob^3+dib^2+{d2-'yodr)bm+{d3-'yiTdr)bm = 2{xi,o+K,o)drbm+i^'ob,i>i+T^ob^o 

(A.39) 

dobt^3 + di6f2 + (^2 - 7o5r)6ci + (c^3 - 7i3^5r)&Co = 2(xi,o + Pb,o)drbco 

+ 2pa,odrbm + W'jj^i + TW>^o , (A.40) 

where we have introduced 

Xi = + i{kWip ■ Vi? + mVip ■ VC) • 

Note that we can write d2 — 7oc?r = — S + |T^Q"(9^ — ik) — (70 — 7o)c^r- 
The solutions of equations (A.38)-(A.40) can be written in the form, 

b^3 = F^aiY),^^^ + G^3(r, dy^^ , G^3 = (A.41) 
b^3 = F^3{Y)e'^^ + G^3{r, d)e'^^ , G^s = (A.42) 
6C3 = F^sirje"'^ + G^spr, ^)e"''> , G^s = . (A.43) 

Projecting equations (A.38)-(A.40) onto e^'^'^ and using (A. 4) gives 

EF^i + 2iaoFi)i = (ps + + \i\io^^ - ~i^9t)F^o - 2(xi,o + J^i,o)F4,o - 2iaJ^Fm 

EFci + W'.F^i = (p3 + ioo3 + izn;"r3 - -f,Tdl)Fco - 2xi,o^co ^ ^Pa^oF^fo " TW^^o , 

since G'^^ = 0, G^;,! = 0, = = Jl'^^^ = 0, Md^o = Jl'^^^ = 0, 7Ig_„ = Jl'g^^ = 
0) o ~ 0- The primes on the G functions denote differentiation with respect 
to r. Note from (A.20) and (A.36), 

2F'" 

In vector form the projected equations are 
LFi = {(p3+ia;3)2/„ + im;"r=^2/„-7iry;'-2xi,„y;}a-ry„ai-2t/;a2, (A.44) 

where L is defined in (A. 25), Fi = (F^i, F^i, F^i)^ and 

ai = [2ia;0;, -^l^n, -WXf , a2 = [-Tli,o^n,\a,o^'o,Pa,o^'of ■ (A.45) 

The primes on indicate derivatives with respect to T. 

We next express the derivatives and terms multiphed by T on the right side 
of (A.44) in terms of paraboUc cy Under functions of different orders by using 
the recurrence relations, 

y'n = {"^yn-l - Vn+l) , = K{ny„-i + Un+l) , (A.46) 
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which arc derived from the parabohc cyhnder function recurrence relations 
(Abramowitz & Stegun 1972). Thus (A. 44) becomes 

3 

LFi = J2' gn,j Vn+j , (A.47) 

where the prime on the summation sign indicates summation over every second 
index. The vectors gn,j are given by 

gn,0 = (P3 + «W3)a 

g„,_3 = in(n - l)(n - 2)(i7on;7n;' - 7^) K'^a 

g„ _i = i[n^7on"7n" + n{n - 2)7^ - inxijK^^a - nnai - nK~^a2 

gn,i = \[{n + l)7on'„7no + {n + 3)7i + ^XiJ^-^a - nai + K'^aa 

gn,3 = i(|7oC/no-7i)«"'a. 

We assume a solution to (A.44) of the form 

3 

i=-3 

By (A.47) the coefficient vectors in (A. 48) are determined from the linear equa- 
tions 

l-n+jh„,j = gn,j , i = 0, ±1, ±3 , (A.49) 

where L„+j is defined in (A.28)(b). When j ^ the determination of b„j is 
straightforward, since J-in+j is invertible, 

1 / Cn+i -2iao 

Kii = 72 « Cn+i I . (A.50) 



n+j el _ t2 1 " 21L, 

We find that 



T — 1 ^ X —1 ^ 

L„±ia = T— a , L„±3a = T^a , 

(A.51) 



and 



^n±l <n \W'J1, ,^^^ + 2moW',n'X,o/^n±l + (^±1 ' en)Pa,oK/in±l 

(A.52) 
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Hence the solution vectors in (A.48) determined from (A. 49) for j ^0 are 

b„,_3 = ^n(n - l)(n - 2)(l7on;7n;' - 7i)^ a (A.53) 

b„,_i = l[n^%IL'^ /U'^ + n{n - 2)j^ - 4nxi,J^a - nL;li(Kai + K-^az) 

(A.54) 

b„,i = -\[{n + l)7on;7n;' + {n + 3)7i + 4xi,o]^a - L-|i(«ai - ^-^a^) 

7o 

(A.55) 

bn,3 = -MtoK/K - 7i)^a. (A.56) 

In the j = case the matrix L„ is singular. It satisfies c^I^n = 0, where 

c = [-n',,^„,of/2n',^rr (A.57) 

and c-'^a = 1. Thus the j = equation in (A.48) furnishes us with the solvability 
condition c^g„,o = 0, which gives us 

P3 + iuj3 = 0, (A.58) 

and which ensures that g„^o = and b„ q is a constant multiple of a. Thus 
the term b„,o Vn can be absorbed into the order solution. Without loss of 
generality we can set b„,o = 0. The particular G solutions can consequently be 
computed, but we omit those details in the interest of space. 

Higher order equations can be reduced to the form above and solved simi- 
larly, though some fortitude is required to weather the algebraic maelstrom that 
ensues. Generally, at some order e^^"^, the projected inhomogeneous equations 
can be manipulated into 

SAT 

L Fjv = ^ iln,j Un+j , 
j=-3N 

which admits a solution of the form 

3JV 
j=-3N 

with the hn,j coefficients determined by linear equations, and the growth rate 
and frequency pn+2 and wjv+2 determined from the singular j = case. 

Appendix B. Numerical evaluation of the asymptotic expressions 

In general the asymptotic approximations (3.24)-(3.27) to the growth rates 
and the eigenfunctions must be evaluated numerically. In this section we de- 
scribe the method of computation. 
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The formula for the growth rate requires the cvahiation of fi, fl', O", W, W 
, W , /Zfc, Mc' /3fc' 0m' Pmk To on the resonant streamhne 4' = '^o- Each 
of these quantities may be evaluated by line integrals along the streamline. 
Moreover a number of their constituent parts (such as Vi? and V^) may also 
be determined by line integrals. The integrals are evaluated numerically using 
the compound trapezoidal rule. The eigenfunctions require the evaluation of 
these quantities, except the /3's, on a (i?, ^')-grid, which must subsequently be 
interpolated onto the (s,z) coordinate system. A simple linear interpolation 
was sufficient. 

Wc obtain Vl by integrating (3.11) and W by averaging. Their derivatives 
may be procured as line integrals using the following technique. The average 
of a function F(\l/,i?) over the curve C$ given by constant ^ can be expressed, 
using (3.11), as 



^ f „ fF 



2tt 



/ V X (^^vm^ • l^rrfrrf6» + niV', (B.l) 



where S^j, is the annular region in the meridional plane bounded by C^j, and a 
smaller fixed ^P-curve Co, which encloses the stagnation point. The quantity K 
is independent of Now by (3.14), rdrdO = {Qr sin 9)~^d^ di!}, and thus 

1 / Vx f:^v,„V-i^dt?rf* + i^. (B.2) 

Differentiation gives 



d^\nj 277 /c^ W J Ursine Q 
where Fi = V x (Fvm/q^) ■ l^/rsin6. Iteration yields 

rf"(F/0)/d*" = K/n , (B.4) 
where F„ is defined inductively for integer n > 1 by 

F„ = V X (F„_iv„,/g2) . l^/rsin0 = -V • ( ^^^f] ■ (B.5) 



(V*)2 

Setting F = 1 gives integral expressions for O', fl" , etc. For azimuthal 
flows of the form W{'d,'i/), F = W gives expressions for W, W , W , etc. 
The integrands rapidly become very complicated with n. The motivation for 
persisting with these complicated expressions is that the numerical integration 
of smooth periodic functions over a period using the compound trapezoidal rule 
is spectrally accurate. 

It soon becomes apparent, however, that a number of the integrands are 
singular in spherical polar coordinates. If Vm = 1^+ Ig, we can write rfi? = 
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s 

Figure B.9: Toroidal coordinates {R, 0). 



flrd9/Vg = Qdr/Vr- The spherical polar components of vanish at two points 
on a VTO-streamline, which can cause singular integrals, specifically those which 
evaluate V?9. To avoid this problem we transform to the toroidal coordinate 
system (i?, O) shown in Figure 9. The point Po(^o,^^o) is the stagnation point 
of the meridional flow. 

Basic trigonometry gives a number of relationships between (r, 9) and {R, 0), 
the most important of which are 

r = ^rl + R'^ + 2roRsm{eo + 6) (B.6) 
^ ^ 1 i?cos(6'o + 6) 

^ = ^o+tan-i \ ° / , B.7 

The {R, 6, (p) system is orthogonal but left-handed. For the flows we examine, 
we take (ro,0o) to be the stagnation point at the centre of the concentric 
curves. The radius ro can be evaluated using the Newton-Raphson method. We 
distinguish between different ^ curves by their largest s- intercept Vs- Each ^ 
curve is thus described by an equation 5'(r, 0) = ^'(r,;, 7r/2). To determine the 
quadrature nodes we divide [0, 27r], corresponding to one period of G, into K 
equal sub-intervals using the nodes O, = 2ni/K, i = 0,1, . . . ,K. To calculate 
the corresponding values of i?, we use the Newton-Raphson method to solve 

^{r{Ri, Qi), e{Ri, Qi)) - *(r„ 7r/2) = 

for Rj given Oi. 

The integrands jJi,, fic, Pk, Prm Pmk and 70 are expressed in terms of: fl, 
W; the R and 9 derivatives of ^, 1?, Vq, Vr; and the ^ and "& derivatives 
of Z; where V„j = VrIr + Vele- The derivatives of ^, and Vu, Vq can 
be found analytically. The quantities dR-d, de^ and d^Z must be determined 
numerically. The angle is given by 

^ = /' dr = r i?* ('-"""'" + f ""Q*^ d9* , (B.8) 
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K 


100 


200 


400 


800 


Vb 


4.63875 


4.63936 


4.63951 


4.63955 


Mc 


-0.20886 


-0.20782 


-0.20756 


-0.20750 


0k 


8.961 


8.964 


8.965 


8.965 


(Vl) 


4.58627 


4.58627 


4.58627 


4.58627 


Prr, (V2) 


4.6719 


4.6718 


4.6718 


4.6718 




0.2256 


0.2278 


0.2284 


0.2285 



Table B.3: The quantities required by the asymptotic theory which converge most slowly with 
K. Here rs = 0.93 and K is the number of subintervals approximating ^fo. 



where an asterisk denotes evaluation upon the — vp* curve. The two sets of 
variables (-R, 8) and (i?*,6*) should not be confused: the asterisked pair are 
dependent on each other while the other pair are independent. With this in 
mind we differentiate the integral in (B.8). Using Leibniz's theorem we obtain 

and 

^ ^ f R* {rosin 9o + R* cos e*)\ de^ R{rosmeo + RcosQ) 

(B.IO) 

Lastly, diiiZ and 9^Zare required for VC- Only the former issues a challenge. 
From the definition of W, 

d^{nZ)= d<^W*M* = d^W*dr -W''&. (B.ll) 
Jo Jo 

Thus from (3.12), 

3.Z = i j\8.W* - "^W^dr + ^{^W- r ) . (B.12) 

This expression is evaluated by converting the ■& integral to an integral over 
and using the formula 

d^W = J-\de^dRW-dRddeW), J = . (B.13) 

Note that as these integrals are not over closed curves, the trapezoidal rule does 
not yield exponential accuracy. The /i's and /3's also require further averaging 
so their convergence is not as fast (as shown in Table 3). 

The asymptotic estimates were computed in MATLAB. The convergence of 
quantities required by the asymptotic theory is shown in Table 3 for different 
numbers K of numerical integration nodes along the chosen streamline, ^o- 
Those quantities that issue from a single integration around the closed stream- 
line converge very rapidly, typically for K = 30. These include ilg and its 4' 
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derivatives, W and its derivatives, and 7q. Thus the quantities flo = 5.3919, 
= 7.2807, n'l = -16.662, Wo = 0.93043, W'^ = 1.4927, W'^' = -2.7786 are 
accurate to the figures shown here for K = 100. However, quantities which are 
evaluated by hne integrals with variable limits converge more slowly, and settle 
down only for K = 800. 
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